feat(rust/sedona-raster-gdal): add RS_Clip#1000
Open
james-willis wants to merge 11 commits into
Open
Conversation
Faithful copy of RS_Clip (src + bench) from the GDAL raster draft PR apache#704, wired into the module list, register, and bench manifest. This commit does NOT compile against current main — it predates several API changes; the reconciliation is isolated into the next commit so the diff there shows exactly what changed from Kristin's original.
Make Kristin's RS_Clip build and pass on today's main: - Add the sedona-proj dependency (CRS reprojection of the clip geometry), plus proj-sys in dev-dependencies for the transform tests. - Drop the removed RasterBandReader; read band bytes via BandRef::nd_buffer().as_contiguous(). - RasterMetadata width/height are now i64 (were u64). - Rewrite the tests onto the RasterSpec harness with self-contained synthetic rasters (no test4.tiff fixture) and RasterStructArray::try_new. Still 2-D only; N-D plane broadcasting follows in the next commit.
A clip is a 2-D (y, x) operation. For an N-D band (extra leading dims such as time), compute the mask and crop window once and apply them to every non-spatial plane, preserving the non-spatial dims and shrinking only the (y, x) extent. Emit the result through the N-D builder (start_raster_nd / start_band_nd); a plain 2-D raster is just the ["y", "x"] case. Adds an end-to-end test over a [time, y, x] raster.
- docs/reference/sql/rs_clip.qmd documenting all five overloads (optional all_touched / no_data_value / crop / lenient) with runnable examples. - Replace the placeholder bench with one that clips a generated raster by a reprojected polygon, over two polygon complexities, via the standard benchmark::scalar harness (needs sedona-testing's criterion feature).
- Output packaging now considers all args, not just [raster, geom]: a per-row band/option column over a scalar raster+geom yields the full N-row array instead of silently collapsing to row 0 (apache#1). - Band 0 (= all bands) reaches clip_raster instead of being clamped to 1 by band.max(1); a negative band errors rather than silently clipping band 1 (apache#2). - Default nodata is the band's own nodata, else the band data type's minimum, never a silent 0.0 that collides with real zero pixels (apache#3). - 'lenient' softens only the no-intersection case; genuine errors (malformed WKB, GDAL failures) always propagate. Dropped the eprintln! (apache#4). - Collapsed the three duplicate nodata/byte-size encoders onto gdal_common::nodata_f64_to_bytes and BandDataType::byte_size() (apache#6). - Dropped the redundant mask zero-init write (GDAL MEM bands are zero-filled) and the double rasterband(1) fetch (apache#8). Adds regression tests for the band-column row-drop, band-0 all-bands, negative-band error, and error-vs-lenient cases.
RS_Clip reads band pixels but wasn't tagged needs_pixels, so the planner never injected RS_EnsureLoaded ahead of it — an OutDb (e.g. Zarr/GeoTIFF-backed) raster would reach the clip with empty band data. Tag needs_pixels (materialize inputs first) and returns_bytes (output is already InDb, don't re-wrap).
…g + real bench - NULL band (or all_touched/crop/lenient) now yields a NULL row instead of silently clipping all bands; no_data_value NULL stays the 'unset' sentinel - reject a no_data_value that can't be represented in the band type instead of saturating it (e.g. -9999 -> 0 on UInt8) - guard the trailing (y, x) dims with is_spatial_dim_pair before clipping - preserve source band names on the clipped output - borrow the band buffer instead of copying it; reserve the plane accumulator - rebuild the benchmark so the raster covers the clip polygon's extent (the clip path is now actually exercised) and drop the copied CRS-tag shim - export rs_clip_udf from lib.rs; document band 0 = all bands and NULL handling
…ll_touched When a strict (lenient=false) clip selects no pixels, the message now depends on all_touched: a plain 'do not intersect' when every touched pixel was already considered (all_touched=true), and a sub-pixel hint pointing at all_touched otherwise. Documents that pixel-center selection can miss sub-pixel geometries.
apply_mask_and_crop copied one pixel at a time via a dynamic-width copy_from_slice, which the compiler can't vectorize. Copy each contiguous crop-window row in one bulk memcpy, then overwrite only the masked-out pixels with nodata (matching the pattern the crop=false path already uses). Add a large-clip benchmark config where this copy dominates.
Contributor
Author
|
@paleolimbot should we hoist the scalar geometry case? Is that a common case we want to optimize? |
Member
|
I think it would be rather difficult to have enough raster tiles that further optimizing the geometry iteration would make a difference given that this is reaching all the way into GDAL and pulling pixels. For what it's worth, I believe the trivial version of this (north up) can be zero copy using views. (Whether that is the best choice for any particular chain of operations is a different story, so maybe not now). |
Contributor
Author
I was thinking this but never finished the views feature |
prantogg
requested changes
Jul 8, 2026
Comment on lines
+601
to
+629
| /// Compute the minimal bounding pixel window that contains all non-zero mask pixels. | ||
| fn compute_crop_window(mask: &[u8], width: usize, height: usize) -> Option<CropWindow> { | ||
| let mut min_col = width; | ||
| let mut max_col = 0usize; | ||
| let mut min_row = height; | ||
| let mut max_row = 0usize; | ||
|
|
||
| for row in 0..height { | ||
| for col in 0..width { | ||
| if mask[row * width + col] != 0 { | ||
| min_col = min_col.min(col); | ||
| max_col = max_col.max(col); | ||
| min_row = min_row.min(row); | ||
| max_row = max_row.max(row); | ||
| } | ||
| } | ||
| } | ||
|
|
||
| if min_col > max_col || min_row > max_row { | ||
| return None; // no non-zero pixels (shouldn't happen if caller checked) | ||
| } | ||
|
|
||
| Some(CropWindow { | ||
| col_off: min_col, | ||
| row_off: min_row, | ||
| width: max_col - min_col + 1, | ||
| height: max_row - min_row + 1, | ||
| }) | ||
| } |
Contributor
There was a problem hiding this comment.
The crop window here is the tight bbox of the selected mask pixels, but the prior art crops to the geometry envelope ∩ raster extent, snapped to the grid — PostGIS ST_Clip (which the module doc cites), gdalwarp -crop_to_cutline, and Sedona Spark's RS_Clip all do the latter, keeping unselected envelope pixels as nodata padding.
Contributor
There was a problem hiding this comment.
There's a performance angle too. Tight crop means the output extent isn't known until after the mask is burned, so every clip pays for a full-raster W×H mask (MEM dataset + the full scans in has_intersection and here) even when the geometry covers a handful of pixels.
Envelope crop knows the window up front, which enables three wins:
- the mask shrinks to the window
- this scan goes away
- and disjoint geometries can be rejected with a cheap envelope check before any GDAL work
Sedona Spark's implementation is structured exactly this way.
…omputed up front The crop window is now the geometry envelope intersected with the raster extent, snapped outward to the pixel grid (PostGIS ST_Clip / gdalwarp -crop_to_cutline semantics), instead of the tight bbox of selected mask pixels. Computing the window before rasterizing rejects disjoint geometries with a cheap envelope check, shrinks the mask to the window, and removes the full-raster crop-window scan.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
3 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Migrates RS_Clip out of the GDAL raster draft (#704) into its own module, with N-D plane-broadcast support, SQL reference docs, and a benchmark.